20. Modelli di equazioni strutturali#
Un modello di equazioni strutturali (SEM) è una tecnica statistica avanzata che combina l’analisi fattoriale confermativa e l’analisi dei percorsi. La modellizzazione delle equazioni strutturali comprende due componenti principali: il modello di misura e il modello strutturale. Il modello di misura può essere rappresentato come un modello fattoriale, mentre il modello strutturale è un modello di analisi dei percorsi che specifica le relazioni tra le variabili latenti in termini di effetti causali, simili all’interpretazione fornita dall’analisi di regressione.
Il processo di analisi SEM è composto dai seguenti passaggi e decisioni:
Costruire un diagramma dei percorsi che mostri il modello di misura e strutturale di interesse.
Identificare il livello di misura per ogni item e verificare le ipotesi distributive.
Assicurarsi che la funzione di adattamento scelta sia basata sui tipi di misura (ad esempio, massima verosimiglianza per misure continue, minimi quadrati ponderato per misure ordinali).
Adattare il modello utilizzando la funzione di adattamento appropriata e valutare l’adattamento del modello utilizzando un insieme di indici.
Una volta stabilito un modello plausibile, interpretare i vari parametri a livello di elemento (ad esempio, saturazioni fattoriali, errori standard, valori R-quadrati, termini di errore, ecc.)
suppressPackageStartupMessages({
library("lavaan")
library("lavaanExtra")
library("lavaanPlot")
library("dplyr")
library("tidyr")
library("knitr")
library("mvnormalTest")
library("semPlot")
})
set.seed(42)
data(HolzingerSwineford1939)
glimpse(HolzingerSwineford1939)
Rows: 301
Columns: 15
$ id <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 11, 12, 13, 14, 15, 16, 17, 18, 19, …
$ sex <int> 1, 2, 2, 1, 2, 2, 1, 2, 2, 2, 1, 1, 2, 2, 1, 2, 2, 1, 2, 2, 1, …
$ ageyr <int> 13, 13, 13, 13, 12, 14, 12, 12, 13, 12, 12, 12, 12, 12, 12, 12,…
$ agemo <int> 1, 7, 1, 2, 2, 1, 1, 2, 0, 5, 2, 11, 7, 8, 6, 1, 11, 5, 8, 3, 1…
$ school <fct> Pasteur, Pasteur, Pasteur, Pasteur, Pasteur, Pasteur, Pasteur, …
$ grade <int> 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, …
$ x1 <dbl> 3.333333, 5.333333, 4.500000, 5.333333, 4.833333, 5.333333, 2.8…
$ x2 <dbl> 7.75, 5.25, 5.25, 7.75, 4.75, 5.00, 6.00, 6.25, 5.75, 5.25, 5.7…
$ x3 <dbl> 0.375, 2.125, 1.875, 3.000, 0.875, 2.250, 1.000, 1.875, 1.500, …
$ x4 <dbl> 2.333333, 1.666667, 1.000000, 2.666667, 2.666667, 1.000000, 3.3…
$ x5 <dbl> 5.75, 3.00, 1.75, 4.50, 4.00, 3.00, 6.00, 4.25, 5.75, 5.00, 3.5…
$ x6 <dbl> 1.2857143, 1.2857143, 0.4285714, 2.4285714, 2.5714286, 0.857142…
$ x7 <dbl> 3.391304, 3.782609, 3.260870, 3.000000, 3.695652, 4.347826, 4.6…
$ x8 <dbl> 5.75, 6.25, 3.90, 5.30, 6.30, 6.65, 6.20, 5.15, 4.65, 4.55, 5.7…
$ x9 <dbl> 6.361111, 7.916667, 4.416667, 4.861111, 5.916667, 7.500000, 4.8…
20.1. Valutazione delle assunzioni#
Un’analisi SEM inizia con una valutazione sulla distribuzione delle variabili endogene. In particolare, dobbiamo valutare l’ipotesi di normalità multivariata. Utilizziamo il pacchetto mvnormalTest per la normalità univariata (test W di Shapiro-Wilk) e multivariata (test di asimmetria e curtosi multivariata di Mardia).
Iniziamo con la normalità univariata.
temp <- HolzingerSwineford1939 |>
select(-c(id, sex, school, ageyr, agemo, grade))
mvnout <- mardia(temp)
## Shapiro-Wilk Univariate normality test
mvnout$uv.shapiro
W p-value UV.Normality
x1 0.9928 0.1582 Yes
x2 0.9697 0 No
x3 0.9523 0 No
x4 0.9827 0.0011 No
x5 0.9769 1e-04 No
x6 0.9538 0 No
x7 0.9908 0.056 Yes
x8 0.9807 4e-04 No
x9 0.9942 0.307 Yes
out = mvnout$mv.test
print(out)
Test Statistic p-value Result
1 Skewness 344.9053 0 NO
2 Kurtosis 3.2344 0.0012 NO
3 MV Normality <NA> <NA> NO
I risultati sia dei test univariati che multivariati indicano che le misure non provengono da distribuzioni univariate o multivariate normalmente distribuite (i risultati ‘No’ nella tabella). Affrontiamo questi problemi nella successiva fase di specifica del modello.
20.2. Specificazione del modello#
La sintassi SEM di lavaan è intuitiva. In primo luogo, definiamo il nostro modello utilizzando la sintassi del modello lavaan e quindi specificare i dettagli tecnici per l’adattamento del modello nella funzione sem.
model <- "
# [-----Latent variables (measurement model)-----]
visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9
# [-----------Mediations (named paths)-----------]
speed ~ visual_speed*visual
textual ~ visual_textual*visual
visual ~ ageyr_visual*ageyr + grade_visual*grade
# [---------Regressions (Direct effects)---------]
speed ~ ageyr + grade
textual ~ ageyr + grade
# [------------------Covariances-----------------]
speed ~~ textual
ageyr ~~ grade
# [--------Mediations (indirect effects)---------]
ageyr_visual_speed := ageyr_visual * visual_speed
ageyr_visual_textual := ageyr_visual * visual_textual
grade_visual_speed := grade_visual * visual_speed
grade_visual_textual := grade_visual * visual_textual
"
Dato che abbiamo riscontrato violazioni dell’ipotesi di normalità multivariata richiesta per i modelli SEM, dobbiamo utilizzare lo stimatore “MLM” come funzione di adattamento. Questo stimatore utilizza una procedura di massima verosimiglianza e fornisce errori standard robusti e una statistica di test scalata di Satorra-Bentler per affrontare i problemi di violazione della normalità multivariata.
È importante notare che i problemi con i dati non normali possono portare a una sottostima degli errori standard, il che può portare a respingere troppo spesso l’ipotesi nulla che un parametro sia zero e ad un’inflazione della statistica chi-quadrato del modello, portando a respingere troppo spesso il modello.
20.3. Bontà di adattamento#
In genere, i ricercatori esaminano le statistiche di adattamento del modello prima di procedere all’interpretazione delle stime dei parametri. L’ipotesi nulla in un’analisi SEM è che la matrice di covarianza implicata o riprodotta dal modello specificato sia statisticamente la stessa della matrice di covarianza di input. Contrariamente al solito test delle ipotesi, speriamo di mantenere l’ipotesi nulla che le due matrici siano statisticamente le stesse.
fit <- sem(model, data=HolzingerSwineford1939, std.lv = TRUE, estimator = "MLM")
Iniziamo a valutare l’adattamento del modello con un test chi-quadrato ottenuto dall’output di lavaan come segue:
out = fitMeasures(fit, c("chisq.scaled", "df.scaled", "pvalue.scaled"))
print(out)
chisq.scaled df.scaled pvalue.scaled
110.247 36.000 0.000
Il rapporto tra chi-quadrato e gradi di libertà è minore di 4, il che è accettabile, anche se indicativo di un fit non eccellente.
La Root Mean Square Error of Approximation (RMSEA) è una misura popolare della discrepanza tra le matrici di correlazione basate sul modello e quelle osservate. Utilizza il chi-quadrato del modello nel suo calcolo ma apporta correzioni in base alla complessità del modello (correzione per la parsimonia) e ha una distribuzione campionaria nota in modo da poter calcolare gli intervalli di confidenza. Abbiamo ottenuto valori RMSEA scalati dall’output di lavaan come segue:
out = fitMeasures(fit, c("rmsea.scaled", "rmsea.ci.lower.scaled", "rmsea.ci.upper.scaled"))
print(out)
rmsea.scaled rmsea.ci.lower.scaled rmsea.ci.upper.scaled
0.083 0.066 0.100
Sono state proposte varie linee guida per l’interpretazione dell’RMSEA: RMSEA <= .05 come soglia per un buon adattamento; RMSEA = .05 - .08 come adattamento ragionevole; RMSEA >= .10 come adattamento scarso. Sulla base della stima puntuale RMSEA ottenuta = .083 e dell’intervallo di confidenza al 90% [.066, .100], concludiamo che il modello ha un adattamento appena accettabile.
Per valutare l’adeguatezza del modello, utilizziamo due ulteriori misure di adattamento: l’Indice di Adattamento Comparativo (CFI) e il Residuo Quadratico Medio Radice Standardizzato (srmr). Il CFI è un indice di adattamento incrementale che confronta il modello considerato con un modello di base ristretto. L’srmr, invece, si basa sulle discrepanze tra le covarianze previste dal modello e le covarianze effettive.
Abbiamo ottenuto i valori scalati del CFI e dell’srmr dall’output di lavaan come segue:
out = fitMeasures(fit, c("cfi.scaled", "srmr"))
print(out)
cfi.scaled srmr
0.925 0.060
Per l’interpretazione di queste misure sono state proposte diverse linee guida. In questo esempio, abbiamo utilizzato come valori soglia CFI >= .90 e srmr <= .08. In base a queste soglie, abbiamo concluso che il valore ottenuto di CFI.scaled = .925 e srmr = .060 fornivano ulteriori prove che il nostro modello si adattava ai dati in modo soddisfacente.
Sulla base di questo insieme di misure di adattamento, possiamo concludere che il modello specificato è plausibile.
20.4. Modello di misurazione#
Dato il modello accettabile, siamo passati all’esame delle varie stime dei parametri. Concentrandoci prima sul modello di misura, abbiamo ottenuto le stime dall’output di lavaan come segue:
standardizedsolution(fit, type = "std.all", se = TRUE, zstat = TRUE, pvalue = TRUE, ci = TRUE) %>%
filter(op == "=~") %>%
select(LV = lhs, Item = rhs, Coefficient = est.std, ci.lower, ci.upper, SE = se, Z = z, "p-value" = pvalue)
| LV | Item | Coefficient | ci.lower | ci.upper | SE | Z | p-value |
|---|---|---|---|---|---|---|---|
| <chr> | <chr> | <dbl> | <dbl> | <dbl> | <dbl> | <dbl> | <dbl> |
| visual | x1 | 0.7775046 | 0.6438504 | 0.9111587 | 0.06819215 | 11.401673 | 0.000000e+00 |
| visual | x2 | 0.4247817 | 0.3137330 | 0.5358305 | 0.05665858 | 7.497218 | 6.528111e-14 |
| visual | x3 | 0.5772006 | 0.4714048 | 0.6829964 | 0.05397843 | 10.693171 | 0.000000e+00 |
| textual | x4 | 0.8556927 | 0.8089949 | 0.9023904 | 0.02382582 | 35.914515 | 0.000000e+00 |
| textual | x5 | 0.8579629 | 0.8172423 | 0.8986834 | 0.02077618 | 41.295515 | 0.000000e+00 |
| textual | x6 | 0.8298582 | 0.7841714 | 0.8755451 | 0.02331005 | 35.600879 | 0.000000e+00 |
| speed | x7 | 0.5988535 | 0.5001953 | 0.6975117 | 0.05033674 | 11.896945 | 0.000000e+00 |
| speed | x8 | 0.7506447 | 0.6645132 | 0.8367763 | 0.04394547 | 17.081277 | 0.000000e+00 |
| speed | x9 | 0.6248603 | 0.5361973 | 0.7135233 | 0.04523704 | 13.813023 | 0.000000e+00 |
Questo output presenta i coefficienti standardizzati (carichi fattoriali) per gli item sulle variabili latenti (LV), gli intervalli di confidenza (ci.lower, ci.upper), gli errori standard (SE), i valori Z (test di Wald) e i valori p che testano l’ipotesi nulla che un coefficiente = 0. I carichi fattoriali variavano da .42 a .86, indicando che l’entità delle relazioni tra gli item e i fattori era adeguata (sebbene non ci siano soglie rigide per i carichi accettabili). Si noti che gli errori standard sono robusti, il che significa che sono corretti per le influenze della non-normalità. Si noti anche che tutti i coefficienti sono statisticamente significativi, il che significa che l’ipotesi nulla che un coefficiente = 0 è respinta.
È anche utile esaminare i valori R2, ovvero i carichi standardizzati al quadrato degli elementi. Nel framework SEM, qualsiasi variabile che ha una freccia rivolta verso di essa è definita come una variabile endogena e avrà un valore R2 associato ad essa. I valori R2 mostrati di seguito per ogni item indicano la percentuale di varianza di quell’item spiegata dalla variabile latente corrispondente. Più alta è la percentuale di varianza di un elemento spiegata dal fattore, migliore è l’elemento nella misurazione del fattore. Abbiamo ottenuto i coefficienti R2 dall’output di lavaan come segue:
parameterEstimates(fit, standardized=TRUE, rsquare = TRUE) %>%
filter(op == "r2") %>%
select(Item=rhs, R2 = est)
| Item | R2 |
|---|---|
| <chr> | <dbl> |
| x1 | 0.6045134 |
| x2 | 0.1804395 |
| x3 | 0.3331605 |
| x4 | 0.7322100 |
| x5 | 0.7361003 |
| x6 | 0.6886647 |
| x7 | 0.3586255 |
| x8 | 0.5634675 |
| x9 | 0.3904504 |
| visual | 0.0913337 |
| textual | 0.3270445 |
| speed | 0.3123506 |
I valori R2 per gli item variano da 0.18 a 0.74 e suggeriscono che gli item abbiano una relazione da piccola a sostanziale con una variabile latente. Si noti che non esiste un limite rigido per i valori R2 accettabili, ma valori oltre .50 sono desiderabili.
20.5. Modello strutturale#
Ora ci concentriamo sul modello strutturale. Otteniamo le stime dall’output di lavaan come segue:
standardizedsolution(fit, type = "std.all", se = TRUE, zstat = TRUE, pvalue = TRUE, ci = TRUE)%>%
filter(op == "~") %>%
select(LV=lhs, Item=rhs, Coefficient=est.std, ci.lower, ci.upper, SE=se, Z=z, 'p-value'=pvalue)
| LV | Item | Coefficient | ci.lower | ci.upper | SE | Z | p-value |
|---|---|---|---|---|---|---|---|
| <chr> | <chr> | <dbl> | <dbl> | <dbl> | <dbl> | <dbl> | <dbl> |
| speed | visual | 0.3795707 | 0.22433706 | 0.53480433 | 0.07920229 | 4.792421 | 1.647810e-06 |
| textual | visual | 0.3660910 | 0.21586372 | 0.51631835 | 0.07664800 | 4.776263 | 1.785824e-06 |
| visual | ageyr | -0.2197484 | -0.37673115 | -0.06276561 | 0.08009472 | -2.743606 | 6.076837e-03 |
| visual | grade | 0.3483088 | 0.19853424 | 0.49808339 | 0.07641700 | 4.558001 | 5.164268e-06 |
| speed | ageyr | 0.1241669 | -0.01806045 | 0.26639422 | 0.07256630 | 1.711082 | 8.706600e-02 |
| speed | grade | 0.2714619 | 0.11145077 | 0.43147309 | 0.08163985 | 3.325116 | 8.838189e-04 |
| textual | ageyr | -0.3852454 | -0.49088915 | -0.27960160 | 0.05390087 | -7.147294 | 8.850698e-13 |
| textual | grade | 0.3229657 | 0.18950793 | 0.45642348 | 0.06809195 | 4.743082 | 2.104911e-06 |
Questo output presenta coefficienti di regressione standardizzati che rappresentano le relazioni tra variabili indipendenti e le variabili dipendenti, intervalli di confidenza (ci.lower, ci.upper), errori standard (SE), valori Z (test di Wald), valori p che testano l’ipotesi nulla che un coefficiente = 0. Un coefficiente di regressione rappresenta la forza della relazione tra una variabile indipendente e il segno rappresenta la direzione della relazione.
20.6. Diagramma di percorso#
nice_lavaanPlot(fit)
lavaanPlot(
model = fit,
node_options = list(shape = "box", fontname = "Helvetica"),
edge_options = list(color = "grey"),
coefs = TRUE
)